CONFORMAL MAPPING METHODS 
FOR INTERFACIAL DYNAMICS^ 



Microstructural evolution is typically beyond the reach of mathematical analysis, but in two 
dimensions certain problems become tractable by complex analysis. Via the analogy be- 
tween the geometry of the plane and the algebra of complex numbers, moving free boundary 
problems may be elegantly formulated in terms of conformal maps. For over half a century, 
conformal mapping has been applied to continuous interfacial dynamics, primarily in models 
of viscous fingering and solidification. Current developments in materials science include 
models of void electro-migration in metals, brittle fracture, and viscous sintering. Recently, 
conformal-map dynamics has also been formulated for stochastic problems, such as diffusion- 
limited aggregation and dielectric breakdown, which has reinvigorated the subject of fractal 
pattern formation. 

Although restricted to relatively simple models, conformal-map dynamics offers unique ad- 
vantages over other numerical methods discussed in this chapter (such as the Level-Set 
Method) and in Chapter 9 (such as the Phase Field Method). By absorbing all geometrical 
complexity into a time-dependent conformal map, it is possible to transform a moving free 
boundary problem to a simple, static domain, such as a circle or square, which obviates the 
need for front tracking. Conformal mapping also allows the exact representation of very 
complicated domains, which are not easily discretized, even by the most sophisticated adap- 
tive meshes. Above all, however, conformal mapping offers analytical insights for otherwise 
intractable problems. 

After reviewing some elementary concepts from complex analysis in §1, we consider the clas- 
sical application of conformal mapping methods to continous-time interfacial free boundary 
problems in §2. This includes cases where the governing field equation is harmonic, bihar- 
monic, or in a more general conformally invariant class. In §3, we discuss the recent use of 
random, iterated conformal maps to describe analogous discrete-time phenonena of fractal 
growth. Although most of our examples involve planar domains, we note in §4 that interfa- 
cial dynamics can also be formulated on curved surfaces in terms of more general conformal 
maps, such as stereographic projections. We conclude in §5 with some open questions and 
an outlook for future research. 



perspective article to appear in the Handbook of Materials Modeling, ed. by S. Yip et al., Vol. I, Ch. 
4, Art. 4.10 (Springer Science and Business Media, 2005). 
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1 Analytic functions and conformal maps 



We begin by reviewing some ba sic concepts from complex analysis found in textbooks suc h 



as 



Churchill and Brownl (|l990| ). For a fresh geometrical perspective, see iNeedhamI (|l997[ ). 



A general function of a complex variable depends on the real and imaginary parts, x and y, 
or, equivalently, on the linear combinations, z = x+iy and z = x—iy. In contrast, an analytic 
function, which is differentiable in some domain, can be written simply as, w = u+iv = f{z). 
The condition, df/dz = 0, is equivalent to the Cauchy-Riemann equations, 

du dv ^ du dv 
dx dy dy dx^ 

which follow from the existence of a unique derivative, 

df du ^ .dv df dv .du 
dx dx dx d{iy) dy dy^ 

whether taken in the real or imaginary direction. 

Geometrically, analytic functions correspond to a special mappings of the complex plane. In 
the vicinity of any point where the derivative is nonzero, f'{z) ^ 0, the mapping is locally 
linear, dw = f'{z) dz. Therefore, an infinitessimal vector, dz, centered at z is transformed 
into another infinitessimal vector, dw, centered aX w = f{z) by a simple complex multi- 
plication. Recalling Euler's formula, (rie*^^)(r2e*^^) = (rir2)e*^^^'^^^\ this means that the 
mapping causes a local stretch by and local rotation by arg/'(z), regardless of the 

orientation of dz. As a result, an analytic function with a nonzero derivative describes a 
conformal mapping of the plane, which preserves the angle between any pair of intersecting 
curves. Intuitively, a conformal mapping smoothly warps one domain into another with no 
local distortion. 

Conformal mapping provides a very convenient representation of free boundary problems. 
The Riemann Mapping Theorem guarantees the existence of a unique conformal mapping 
between any two simply connected domains, but the challenge is to derive its dynamics 
for a given problem. The only constraint is that the conformal mapping be univalent, or 
one-to-one, so that physical fields remain single-valued in the evolving domain. 
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2 Continuous interfacial dynamics 



2.1 Harmonic fields 



Most applications of conformal mapping involve harmonic funtions, which are solutions to 
Laplace's equation, 

= 0. (3) 

From Eq. it is easy to show that the real and imaginary parts of an analytic function 
are harmonic, but the converse is also true: Every harmonic is the real part of an analytic 
function, = Re the complex potential. 

This connection easily produces new solutions to Laplace's equation in different geometries. 
Suppose that we know the solution, = Re<I>(w), in a simply connected domain in 

the w-plane, Q^, which can be reached by conformal mapping, w = f{z,t), from another, 
possibly time-depedent domain in the z-plane, ^^(t). A solution in Qz{t) is then given by 

(j){z, t) = Re = Re ^{f{z, t)) (4) 

because $(/(z)) is also analytic, with a harmonic real part. The only caveat is that the 
boundary conditions be invariant under the mapping, which holds for Dirichlet ((/)=constant) 
or Neumann (n ■ = 0) conditions. Most other boundary conditions invalidate Eq. Q and 
thus complicate the analysis. 

The complex potential is also convenient for calculating the gradient of a harmonic function. 
Using Eqs. ((T)) and 0, we have 

V.^ = ^ + i|^ = F, (5) 

OX oy 

where Vz is the complex gradient operator, representing the vector gradient, V, in the 
2;-plane. 

Viscous fingering and solidification: The classical application of conformal-map dynam- 
ics is to Laplacian growth, where a free boundary, Bz{t), moves with a (normal) velocity, 

dz 

v = ^ocV0, (6) 
proportional to the gradient of a harmonic function, (f), which vanishes on the boundary 



( HowisonI ()l992n ). Conformal mapping for Laplacian growth was introduced independently 



by Polubarinova-Kochina and Galin in 1945 in the context of groundwater flow, where 
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Figure 1: Exact solutions for Laplacian growth, a simple model of viscous fingering: (a) 
a Saff man- Taylor finger translating down an infinite channel, showing iso-pressure curves 
(dashed) and streamlines (solid) in the viscous fluid, and (b) the evolution of a perturbed 
circular bubble leading to cusp singularities in finite time. [Courtesy of Jaehyuk Choi.] 



is the pressure field and u = {k/ri)'V(f) is the velocity of the fluid of viscosity, t], in a 
porous medium of permeability, k, according to Darcy's law. Laplace's equation follows 
from incompressibility, V ■ u = 0. The free boundary represents an interface with a less 
viscous, immiscible fluid at constant pressure, which is being forced into the more viscous 
fluid. 

In physics, Laplacian growth is viewed as a fundamental model for pattern formation. It also 
describes viscous fingering in Hele-Shaw cells, where a bubble of fluid, such as air, displaces 
a more viscous fluid, such as oil, in the narrow gap between parallel flat plates. In that case, 
the depth averaged velocity satisfies Darcy's law in two dimensions. Laplacian growth also 
describes de ndritic solidification in t he limit of low undercooling, where (p is the temperature 



in the melt (jCummings et al.l (|l999l )). 



To illustrate the derivation of conformal-map dynamics, let us consider viscous fingering in a 
channel with inpenetrable walls, as shown in Fig. ^a). The viscous fluid domain, Qz{t), lies 
in a periodic horizontal strip, to the right of the free boundary, Bz{t), where uniform flow of 
velocity, U, is assumed far ahead of the interface. It is convenient to solve for the conformal 
map, z = g{w,t), to this domain from a half strip, Rew > 0, where the pressure is simply 
linear, = ReUw/fi. We also switch to dimensionless variables, where length is scaled to a 
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characteristic size of the initial condition, L, pressure to UL/fi, and time to L/U. 

Since = 1 in the half strip, the pressure gradient at a point, z = g{w, t), on the physical 
interface is easily obtained from Eq. (|Hn|l : 



V.^ = ^=(|^l (7) 



dz \dw 

where w = f{z,t) is the inverse mapping (which exists as long as the mapping remains 
univalent). Now consider a Lagrangian marker, z{t), on the interface, whose pre-image, 
w{t), lies on the imaginary axis in the w-plane. Using the chain rule and Eq. (jZj), the 
kinematic condition, Eq. (jHl), becomes. 



dz dg dg dw f dg 
dt dt dw dt V dw 



Multiplying by dg/dw ^ 0,this becomes 



dg dg ^ 
dw dt 



dg 



dw 



Since the pre-image moves along the imaginary axis, lie [dw/dt) = 0, we arrive at the 
Poluharinova-Galin equation for the conformal map: 



Ref-^^1=1, forRew = 0. (10) 
\dw dt) ' ^ ' 

From the solution to Eq. ()10|1 . the pressure is given hy cj) = R.ef{z,t). Note that the 
interfacial dynamics is nonlinear, even though the quasi-steady equations for are linear. 

The best-known solutions are the Saff man- Taylor fingers, 

g{w, t) = j+w + 2{l-X) log(l + e-'") (11) 



which translate at a constant velocity, A ^, without changing their shape (jSaffman and Tavloi 



( 1958^ ) . Note that ()11|) is a solution to the fingering problem for all choices of the parameter 
A. This parameter specifies the finger width and can be chosen arbitrarily in the solution 
(fTT|l . In experiments however, it is found that the viscous fingers that form are well fit by 
a Saff man- Taylor finger filling precisely half of the channel, that is with A = 1/2, as shown 
in Fig. ^a). Why this happens is a basic problem in pattern selection, which has been the 
focus of much debate in the literature over the last 25 years. 
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Figure 2: Numerical simulation of viscous fingering, starting from a three-fold perturbation 
of a circular bubble. The only difference with the Laplacian-growth dynamics in Fig. ^b) is 
the inclusion of surface tension, which prevents the formation of cusp singularities. [Courtesy 
of Michael Siegel.] 



To understand this problem, note that the viscous finger solutions (fTTj) do not include any of 
the effects of surface tension on the interface between the two fluids. The intriguing pattern 
selection of the A = 1/2 finger has been attributed to a singular perturbation effect of small 
surface tension. Surface tension, 7, is a significant complication because it is described by a 
non-conformally-invariant boundary condition. 



7K, for z G Bz{t) 



(12) 



where k is the local interfacial curvature, entering via the Young-Laplace pressure. Small 
surface tension c an be treated analytically as a singular per turbation to gain insights into 
pattern selection ( Kruskal and Segui (1991); Tanveer (|l993a )). Since surface tension effects 
are only significant at points of high curvature k in the interface, and given that the finger in 
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Fig. ^a) is very smooth with no such points of high curvature, it is surprising that surface 
tension acts to select t he finger width. Indeed, the viscous fingering problem has been shown 
to be full of surprises (ITanveeii (j2000|)). 



In a radial geometry, the univalent mapping is from the exterior of the unit circle, |w 
the ex terior of a finite bubble penetrating an infinite viscous liquid. iBensimon and Shraiman 
f)l984 ) introduced a pole dynamics formulation, where the map is expressed in terms of its 
zeros and poles, which must lie inside the unit circle to preserve univalency. They showed 
that Laplacian growth in this geometry is ill-posed, in the sense that cusp-like singularites 
occur in finite time (as a ze ro hits the unit circle) for a broad class of initial conditions, as 
illustrated in Fig. mb). fSee lHowisonI (|l992| ) for a simple, general proof due to Hohlov.) This 
initiated a large body of work on how Laplacian growth is "regularized" by surface tension 
or other effects in real systems. 

Despite the analytical complications introduced by surface tens i on, several exa ct steady 
solutions with non-zero surface tension are known (jKadanofa ()l99fll ): ICrowdvl ()2nnfl( )). Surface 
tension can also be incorporated into nurnerical simulations based on the same conformal- 
mapping formalism (IMaclean and SaffmanI (|l98l|)), which show how cusps are avoided by the 
formation of new fingers ([Pai et al' (1991)). For example, consider a three- fold perturbation 
of a circular bubble, whose exact dynamics without surface tension is shown in Fig. ^b). 
With surface tension included, the evolution is very similar until the cusps begin to form, 
at which point the tips bulge outward and split into new finger s, as shown in Fig. El Th is 
process repeats itself to produce a complicated fractal pattern (|Bunde and Havlin 
which curiously ressembles the diffusion-limited particle aggregates discussed below in §3. 

Density-driven instabilities in fluids: An important class of problems in fluid mechanics 
involves the nonlinear dynamics of an interface between two immiscible fluids of different 
densities. In the presence of gravity, there are some familiar cases. Deep-water waves involve 
finite disturbances (such as steady "Stokes waves'") in the interface between lighter fluid 
(air) over a heavier fluid (water). With an inverted density gradient, the Rayleigh- Taylor 
instability develops when a heavi er fluid lies abov e a lighter fluid, leading to large plumes of 
the former sinking into the latter. iTanveeii (|l993bf ) has used conformal mapping to analye the 
Rayleigh- Taylor instability and has provided evidence to associate the formation of plumes 
with the approach of various conformal mapping singularities to the unit circle. 

A related problem is the Richtmyer-Meshkov instability, which occurs when a shock wave 
passes through an interface between fluids of different densities. Behind the shock, narrow 
fingers of the heavier fluid penetrate into the lighter fluid. The shock wave usually passes so 
quickly that compressibility only affects the onset of the instability, while the nonlinear evo- 
lution occurs much faster than the development of viscous effects. Therefore, it is reasonable 
to assume a potential flow in each fluid region, with randomly perturbed initial velocities. 
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Figure 3: Conforma l-map dynamics for the stro ngly nonlinear regime of the Richtmyer- 
Meshkov instability ( Yoshikawa and Balk ( 1999f )). [Courtesy of Toshio Yoshikawa and 
Alexander Balk.l 



Although real plumes roll up in three dimensions and bubbles can form, conformal mapping 
in two dimensions still provides some insights, with direct relevance for shock tubes of high 
aspect ratio. 

A simple conformal- mapping analysis is possible for the case of a large density contrast, where 
the lighter fluid is assumed to be at uniform pressure. The Richtmyer-Meshkov instability 
(zero-gravity limit) is then similar to the Saffman- Taylor instability, except that the total 
volume of each fluid is fixed. A periodic interface in the y direction, analogous to the channel 
geometry in Fig.^ can be described by the univalent mapping, z = g{w, t), from the interior 
of the unit circle in the mathematical w plane to the interior of the heavy-fluid finger in the 
physical z plane. 

Zakharov (|l968) introduced a Hamiltonian formulation of the interfacial dynamics in terms 



of this conformal map, taking into account kinetic and potential energy, but not surface 
tension. One way to derive equations of motion is to expand the map in a Taylor series, 

oo 

g{w,t) = log w + '^an{t)w"', |w| < 1. (13) 

n=0 
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(The logw term first maps the disk to a periodic half strip.) On the unit circle, w = e*^, the 
pre-image of the interface, this is simply a complex Fourier series. The Taylor coefficients, 
a„(t), act generalized coordinates describing n-fold shape perturbations within each period, 
and their time derivatives, an{t), act as velocities or momenta. Unfortunately, truncating 
the Taylor series results in a poor description of strongly nonlinear dynamics because the 
conformal map begins to vary wildly near the unit circle. 



An alternate approach used by lYoshikawa and Balkl ()l999l ) is to expand in terms ressembling 
Saff man- Taylor fingers. 



N 



g{w, t) = \ogw + b{t) - bn{t) log(l - Xn{t)w) (14) 

n=l 

which can be viewed as a resummation of the Taylor series in Eq. (|13p. As shown in Fig. |21 
exact solutions exist with only a finite number terms in the finger expansion, as long as 
the new generalized coordinates, Xnit), stay inside the unit disk, |A„| < 1. This example 
illustrates the importance of the choice of shape functions in the expansion of the conformal 
map, e.g. w"' versus log(l — Xnw). 

Void electro-migration in metals: Micron-scale interconnects in modern integrated cir- 
cuits, typically made of aluminum, sustain enormous currents and high temperatures. The 
intense electron wind drives solid-state mass diffusion, especially along dislocations and grain 
boundaries, where voids also nucleate and grow. In the narrowest and most fragile inter- 
connects, grain boundaries are often well separated enough that isolated voids migrate in a 
fairly homogeneous environment due to surface diffusion, driven by the electron wind. Voids 
tend to deform into slits, which propagate across the interconnect, causing it to sever. A 
theory of void electro-migration is thus important for predicting reliability. 



In the simplest two-dimensional model (|Wang et al.l (|l996| )). a cylindrical void is modeled as 



a deformable, insulating inclusion in a conducting matrix. Outside the void, the electrostatic 
potential, (p, satisfies Laplace's equation, which invites the use of conformal mapping. The 
electric field, E = — V0, is taken to be uniform far away and wraps around the void surface, 
due to a Neuman boundary condition, n • E = 0. 

The difference with Laplacian growth lies in the kinematic condition, which is considerably 
more complicated. In place of Eq. (jU)), the normal velocity of the void surface is given by 
the surface divergence of the surface current, j, which takes the dimensionless form, 

n-v = — = ^ + ^ C15) 
ds ds"^ ds"^ 

where s is the local arc-length coordinate and x is a dimensionless parameter comparing 
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(a) 



(b) 



Figure 4: Numerical c onformal-mapping simulations of the electromigration of voids in 
aluminum interconnects ( Wang et al. ( 1996[ l). (a) A small shape perturbation of a cylindrical 
void decaying (above) or deforming into a curved slit (below), depending on a dimensionless 
group, comparing the electron wind to surface-tension gradients, (b) A void evolving with 
anisotropic surface diffusivity (x = 100, = 100, m = 3). [Courtesy of Zhigang Suo.] 



surface currents due to the electron wind force (first term) and due to gradients in surface 
tension (second term). This moving free boundary problem somewhat resembles the viscous 
fingering problem with surface tensio n, and it admits analogous finger solutions, albeit of 
width 2/3, not 1/2 (|Ben Amail ^99^). 



To describe the evolution of a singly connected void, we consider the conformal map, z = 
g{w,t), from the exterior of the unit circle to the exterior of the void. As long as the map 
remains univalent, it has a Laurent series of the form. 



g{w, t) = Ai{t)w + Ao{t) + A_n{t)w-'', for\w\ > 1, 



(16) 



n=l 



where the Lauren t coefficients, A„ (t), are now the generalized coordinates. As in the case of 
viscous fingering (jHowisonl ()l992i )). a hierarchy of nonlinear ordinary di fferential equations 
(ODEs) for these coordinates can be derived. For void electromigration. IWang et al.l (|l996r i 
start from a variational principle accounting for surface tension and surface diffusion, using 
a Galerkin procedure. They truncate the expansion after 17 coefficients, so their numerical 
method breaks down if the void deforms significantly, e.g. into curved slit. Nevertheless, as 
shown in Fig,. ID^a), the numerical method is able to capture essential features of the early 
stages of strongly nonlinear dynamics. 
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In the same regime, it is also possible to incorporate anisotropic surface tension or surface 
mobility. The latter involves multiplying the surface current by a factor (1 + cos ma), 
where a is the surface orientation in the physical z plane, given aX z = (7(e*^,t), by 



dg 



(17) 



Some results are shown in Fig. ^b), where the void develops dynamical facets. 



Quadrature domains: We end this section by commenting on some of the mathematics 
underlying the existence of exact solutions to continuous-time Laplacian-growth problems. 
Significantly, much of this mathematics carries over to problems in which the governing field 
equation is not necessarily harmonic, as will be seen in the following section. 

The steadily-translating finger solution (jllj) of Saffman and Taylor turns out to be but one of 
an increasingly large number of known exact solutions to the standard Hele-Shaw problem. 
SaffmanI ()l959| ) him self identi fi ed a c lass of unsteady finger-like solutions. This solution was 



later generalized by iHowisonI ()l986| ) to solutions involving multiple fingers exhibiting such 
phenomena as tip-splitting where a single finger splits into two (or more) fingers. It is even 
possible to find exact solu tions to the more real i stic ca se where there is a second interface 
further down the channel (jCrowdy and Tanveen ()2n04l )) which must always be the case in 
any experiment. 

Besides finger-like solutions which are characterized by time-evolving conformal mappings 
having logarithmic branch-point singularities, other exact solutions, where the conformal 
mappings are rational functions with time- evolving pol e s and zeros, were first identified by 
Polubarinova-Kochina and Galin in 1945. iRichardsonI ()l98lh later rediscovered the latter 
solutions while simultaneously presenting important new theoretical connections between 
the Hele-Shaw problem and a class of planar domains known as quadrature domains. The 
simplest example of a quadrature domain is a circular disc D of radius r centred at the origin 
which satisfies the identity 




h(z)dxdy = 7rr^/i(0) 



D 



where h{z) is any function analytic in the disc (and integrable over it). Equation (jlSp . which 
is known as a quadrature identity since it holds for any analytic fun ction h(z), is simply a 
statement of the well-known mean-value theorem of complex analysis (jCarrier et al.l (|l966| )). 
A more general domain D, satisfying a generalized quadrature identity of the form 




h{z)dxdy 



D 



N nfe-1 

5^5^c,,/^(^-)(z,) 

k=l j=0 



(19) 
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is known as a quadrature domain. Here, {zk G C} is a set of points inside D and h^^\z) 
denotes the j-th derivative of h{z). If one makes the choice h{z) = in (|19p the resulting 
integral quantities have become known as the Richardson moments of the domain. Richard- 
son showed that the dynamics of the Hele-Shaw problem is such as to preserve quadrature 
domains. That is, if the initial fluid domain in a Hele-Shaw cell is a quadrature domain at 
time t = 0, it remains a quadrature domain at later times (so long as the solution does not 
break down). This result is highly significant and provides a link with many other areas of 
mathematics including potential t heory, the notion o f balayage, algebraic curves, Schwarz 
functi ons and Cauchy transforms. Richardson discusses many of these connections 

while IVarchenko and Etingoa (|l992i l provide a more g enera l overview of the various mathe- 
matical implications of Richardson's result. IShapird ()l992n gives more general background 
on quadrature domain theory. 



It is a well-known result in the theory of quadrature domains (|Shapird (|l992| )) that simply- 
connected quadrature domains can be parametrized by rational function conformal mappings 
from a unit circle. Given Richardson's result on the preservation of quadrature domains, 
this explains why Polubarinova-Kochina and Galin were able to find time-evolving rational 
function conformal mapping solutions to the Hele- Shaw problem. It also underlies the pole 
dynamics results of lBensimon and ShraimanI ()l984^ . But Richardson's result is not restricted 
to simply-connected domains; multiply- connected quadrature domains are also preserved by 
the dynamics. Physically this corresponds to time-evolving fluid domains containing mul- 
tiple bubbles of air. Indeed, motivated by such matters, recent research has focused on 
the analyti cal construct i on of r nultiply-connected quadratur e domains using conformal map- 
ping ideas (|Richardsonl (|2001l l. lCrowdv and Marshalll (|2004l )). In the higher-connected case, 
the conformal mappings are no longer simply rational functions but are given by conformal 
maps that are automorphic functions (or, meromorphic functions on compact Riemann sur- 
faces). The important point here is that understanding the physical problem from the more 
general perspective of quadrature domain theory has led the way to the unveiling of more 
sophisticated classes of exact conformal mapping solutions. 



2.2 Bi-Harmonic fields 



Although not as well known as conformal mapping of harmonic functions, there is also a 
substantial literature on complex-variable methods to solve the bi-harmonic equation, 



0, 



(20) 



which arises in two-dimensional elasticity (|Muskhelishvilil(|l953[ )) and fluid mechanics (iBatcheloij 
( 19671 )). Unlike harmonic functions, which can be expressed in terms of a single analytic 
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Figure 5: Evolution of the solution of iHopperl ()l99fl| ) for the coalescence of two equal blobs 
of fluid under the effects of surface tension. 



function (the complex potential), bi-harmonic func tions can be expressed in terms of two 



3xpr 

analytic functions, f{z) and g{z), in G ours at form (jCarrier et al.l ()l96fi) ): 



^P{z,z) = lm[zf{z)+giz)] (21) 

Note that ip is no longer just the imaginary part of an analytic function g{z) but also 
contains the imaginary part of the non-analytic component zf{z). A difficulty with bi- 
harmonic problems is that the typical boundary conditions (see below) are not conformally 
invariant, so conformal mapping does not usually generate new solutions by simply a change 
of variables, as in Eq. (@)). Nevertheless, the Goursat form of the solution, Eq. (|2ip. is a 
major simplification, which enables analytical progress. 

Viscous sintering: Sintering describes a process by which a granular compact of particles 
(e.g. metal or glass) is raised to a sufficiently large temperature that the individual particles 
become mobile and release surface energy in such a way as to produce inter-particulate bonds. 
At the start of a sinter process, any two particles which are initially touching develop a thin 
"neck" which, as time evolves, grows in size to form a more developed bond. In compacts in 
which the packing is such that particles have more than one touching neighbour, as the necks 
grow in size, the sinter body densifies and any enclosed pores between particles tend to close 
up. The macroscopic material properties of the compact at the end of the sinter process 
depend heavily on the degree of densification. In industrial application, it is crucial to be 
able to obtain accurate and reliable estimates of the time taken for pores to close (or reduce 
to a sufficiently small size) within any given initial sinter body in order that industrial sinter 
times are optimized without compromising the macroscopic properties of the final densified 
sinter body. 

The fluid is modelled as a region D[t) of very viscous, incompressible fluid, in which the 
velocity fleld, 

n={u,v) = {iljy,-i:,). (22) 
is given by the curl of an out-of-plane vecto r, whose ni a gnitu de is a stream function, ipi^x, y, t). 



which satisfles the bi-harmonic equation (jBatcheloii ()l967l )). On the boundary dD{t), the 
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tangential stress must vanish and the normal stress must be balanced by the uniform surface 
tension effect, i.e., 

— pui + 2/ieij = Turii, (23) 

where p is the fluid pressure, /i is the viscosity, T is the surface tension parameter, k is the 
boundary curvature, rij denotes components of the outward normal n to dD{t) and Cij is 
the usual fluid rate-of-strain tensor. The boundary is time-evolved in a quasi-steady fashion 
with a normal velocity, Vn, determined by the same kinematic condition, = u ■ n, as in 
viscous fingering. 

In terms of the Goursat functions in (j2H) - which are now generally time-evolving - the 
stress condition takes the form 

fiz, t) + zj{z, t) + t) = --Z, (24) 

where again s denotes arclength. Once f{z,t) has been determined from (j23)) the kinematic 
condition ^ 

lm[ztZs] = Im[-2/(z, t)zs] - - (25) 

is used to time-advance the interface. 



A significant contribution was made bv iHopperl ()l990| ) who showed, using complex variable 



methods based on the decomposition ()21|) . that the problem for the surface-tension driven 
coalescence of two equal circular blobs of viscous fluid can be reduced to the evolution of a 
rational function conformal map, from a unit w-circle, of the form 

R(t)w , , 

9{w,t)= ^ ^\ 26 

The two time-evolving parameters R{t) and a{t) satisfy two coupled nonlinear ODEs. Figure 
shows a sequence of shapes of the two coalescing blobs computed using Hopper's solution. 
At large times, the configuration equilibrates to a single circular blob. 

While Hopper's coalescence problem provides insight into the growth of the interparticle 
neck region, there are no por es in this confi guration and it is natural to ask whether more 



general exact solutions exist. ICrowdvl (|1999I ) reappraised the viscous sintering problem and 
showed, in direct analogy with Richardson's result on Hele-Shaw flows, that the dynamics 
of the sintering problem is also such as to preserve quadrature domains. As in the Hele- 
Shaw problem, this perspective paved the way for the identification of new exact solutions, 
generalizing ()2fj|l . for the evolution of doubly-connected fluid regions. Figure El shows the 
shrinkage of a pore enclosed by a typical "unit" in a doubly-connected square packing of 
touching near-circular blobs of viscous fluid. This calculation employs a conformal mapping 
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Figure 6: The coalescence of fluid blobs and collapse of cylindrical pore s in a ni o del o f 
viscous sintering. This sequence of images shows an analytical solution by Crowdv ( 2003[ l 
using complex-variable methods. 



to the doubly- connected fluid regi on (which is no longer a rational function but a more general 
automorphic function) derived bv lCrowdvl (|2003h and, in the same spirit as Hopper's solution 
(j2Sl); requires only the integration of three coupled nonlinear ODEs. The fluid regions 
shown in Figure El are all doubly-connected quadrature domains. iRichardsonI ()2000[ ) has also 
considered similar Stokes flow problems using a different conformal mapping approach. 

Pores in elastic solids: Solid ela sticity in two dinaensio ns is also governed by a bi-harmonic 
function, the Airy stress function (iMuskhelishvilil (|l953|)). Therefore, the stress tensor, aij, 
and the displacement fleld, Ui, may be expressed in terms of two analytic functions, f{z) 
and g{z): 



0"22 - CTll 



+ i(Ti2 



Y 



1 + 



Ui + iU2) 



f\z) + f'{z) 
zfiz) + g'{z) 
i^f{z) -zf'{z) -g{z) 



(27) 
(28) 
(29) 



where Y is Young's modulus, u is Poisson's ratio, and k = (3 — z/)/(l + z/) for plane stress and 
K = 3 — 4i> for plane strain. As with bubbles in viscous flow, the use of Goursat functions 
allows conformal mapping to be applied to bi-harmonic free boundary problems in elastic 
solids, without solving explicitly for bulk stresses and strains. 



For example, IWang and Sud (|l997f l have simulated the dynamics of a singly-connected pore 
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by surface diffusion in an infinite stressed elastic solid. As in the case of void electromigration 
described above, they solve nonlinear ODEs for the Laurent coefficients of the conformal map 
from the exterior of the unit disk, Eq. Under uniaxial tension, there is a competition 

between surface tension, which prefers a circular shape, and the applied stress, which drives 
elongation and eventually fracture in the transverse direction. The numerical method, based 
on the truncated Laurent expansion, is able to capture the transition from stable elliptical 
shapes at small applied stress to the unstable growth of transverse protrusions at large 
applied stress, although naturally it breaks down when cusps ressembling crack tips begin 
to form. 



2.3 Non-harmonic conformally invariant fields 

The vast majority of applications of conformal mapping fall into one of the two classes above, 
involving harmonic or bi-harmonic functions, where the connections with analytic functions, 
Eqs. dH) and (PT|). are cleverly exploited. It turns out, however, that conformal mapping can 
be applied ius t as easi l y to a broad class of problems involving non-harmonic fields, recently 



discovered bv iBazand (j200J). Of course, in planar geometry, the conformal map itself is 
described by an analytic function, but the fields need not be, as long as they transform in a 
simple way under conformal mapping. 

The most convenient fields satisfy conformally invariant partial differential equations (PDEs), 
whose forms are unaltered by a conformal change of variables. It is straightforward to trans- 
form PDEs under a conformal mapping of the plane, w = f{z),hj expressing them in terms 
of complex gradient operator introduced above, 

d d d 

^^ = 7r + *ir = 2^' (30) 

ox ay oz 

which we have related to the 1 partial derivative using the Cauchy-Riemann equations, 
Eq. (^. In this form, it is clear that Vz/ = if and only if f{z) is analytic, in which case 
V zf = 2/'. Using the chain rule, also obtain the transformation rule for the gradient, 

V, = fV^ (31) 

To apply this formalism, we write Laplace's equation in the form, 

= Re V, = V, V,0 = 0, (32) 

which assumes that mixed partial derivatives can be taken in either order. (Note that 
a ■ b = Rea6.) The conformal invariance of Laplace's equation, V^V^0 = 0, then follows 
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from a simple calculation, 



V,V, = iVJ')\/^ + l/f V^V^ = \ff V^V» (33) 

where Vz/' = because /' is also analytic. As a result of conformal invariance, any harmonic 
function in the w plane, (piw), remains harmonic in the z plane, (j){f{z)), after the simple 
substitution, w = f{z). We came to the same conclusion above in Eq. (jH), using the 
connection between harmonic and analytic functions, but the argument here is more general 
and also applies to other PDEs. 

The bi-harmonic equation is not conformally invariant, but some other equations - and 
systems of equations - are. The key observation is that any "product of two gradients" 
transforms in the same way under conformal mapping, not only the Laplacian, V ■ V0, but 
also the term, V0i ■ V</)2 = Re (V0i) 5 which involves two real functions, 01 and (f)2- 

Re (Vz0i) V,02 = I/? Re (V^0i) V^02. (34) 

(Todd Squires has since noted that the operator, V0i x V02 = Im (V(/)i) V02, also transforms 
in the same way.) These observations imply the conformal invariance of a broad class of 
systems of nonlinear PDEs: 

N / N N \ 

i=l \ j=i j=i+l J 

where the coefficients ai{4>), aij{4>), and bij{4>) may be nonlinear functions of the unknowns, 
(f) = (01, 02, . . . , (Pn), but not of the independent variables or any derivatives of the unknowns. 

The general solutions to these equations are not harmonic and thus depend on both z and 
1. Nevertheless, conformal mapping works in precisely the same way: A solution, (f){w,w), 
can be mapped to another solution, (f){f{z), f{z)), by a simple substitution, w = f{z). This 
allows the conformal mapping techniques above (and below) to be extended to new kinds of 
moving free boundary problems. 

Transport-limited growth phenomena: For physical applications, the conformally in- 
variant class, Eq. (jH^j) . includes the important set of steady conservation laws for gradient- 
driven flux densities, 

(9c • 

^ = V-Fi = 0, Fi = CiUi-Di{ci)Vci, u,, cxV0 (36) 
ot 

where {q} are scalar fields, such as chemical concentrations or temperature, {Di{ci)} are 
nonlinear diffusivities, {uj} are irrotational vector fields causing advection, and is a poten- 
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Figure 7: The exact self-similar solution, Eq. ()40j) . for continuous advection-diffusion-limited 
growth in a uniform background potential flow (yellow streamlines) at the dynamical fixed 
point (Pe = oo). The concentration field (color contour plot) is shown for Pe = 100. 
[Courtesy of Jaehyuk Choi.] 



tial (jBazantl 1)20041 )). Physical examples include advect ion- diffusion, where is the harmonic 
velocity potential, and electrochemical transport, where is the non-harmonic electrostatic 
potential, determined implicitly by electro-neutrality. 

By modifying the classical methods described above for Laplacian growth, conformal-map 
dynamics can thus b e formulated for more general, transport-limited growth phenomena 
fiBazant et al.l ()2003| )). The natural generalization of the kinematic condition, Eq. (jHl), is 
that the free boundary moves in proportion to one of the gradient-driven fluxes with velocity, 
V oc Fi. For the growth of a finite filament, driven by prescribed fluxes and/or concentrations 
at infinity, one obtains a generalization of the Polubarinova-Galin equation for the conformal 
map, z = g{w, t), from the exterior of the unit disk to the exterior of growing object. 



Re{wg' gt) = (j{w,t), on |w| = 1 



(37) 



where a{w,t) is the non-constant, time-dependent normal flux, n ■ Fi, on the unit circle in 
the mathematical plane. 

Solidification in a fluid flow: A special case of the conformally invariant equations ()H3|l 
has been known for almost a century: steady advection-diffusion of a scalar field, c, in a 
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potential flow, u. The dimensionless PDEs are 



PeS^=(^ + S) (39) 



Pe u ■ Vc = V^c, u = V0, = (38) 

where we have introduced the Peclet number, Pe = UL/D, in terms of a characteristic 
length, L, velocity, f/, and diffusivity, D. In 1905, Boussinesq showed that Equation 
takes a simpler form in streamline coordinates, {(f),ip), where $ = + i-^ is the complex 
velocity potential: 

dc ( (9^c d'^c 

because advection (the left hand side) is directed only along streamlines, while diffusion 
(the right hand side) also occurs in the transverse direction, along isopotential lines. From 
the general perspective above, we recognize this as the conformal mapping of an invariant 
system of PDEs of the form (jHBj) to the complex $ plane, where the flow is uniform and any 
obstacles in the flow are mapped to horizontal slits. 

Streamline coordinates form the basis for Maksimov's method for interfacial growth by 
steady advection-diffusion in a background potential flow, w hich has been applied to freez- 



i ng in g roundwater flow and vap or deposition on textile flbers (jKornev and MukhamaduUina 



( 1994 ): ICummings et al. 1 (|l999)). The growing interface is a streamline held at a flxed con- 



centration (or temperature) relative to the flowing bulk fluid at inflnity. This is arguably 
the simplest growth model with two competing transport processes, and yet open questions 
remain about the nonlinear dynamics, even without surface tension. 

The normal flux distribution to a flnite absorber in a unifor r n backg round flow, a{w,t) in 



Eq. (j37|) is well known, but rather complicated (jChoi et al.l (|2004b| )). so it is replaced by 



asymptotic approximations for analytical work, such as a ~ 2A/Pe /7rsinfg/2) as Pe 



which is the flxed point of the dynamics. In this important limit, IChoi et al. 
found an exact similarity solution. 



(200M) have 



g{w,t) = A,{t)^/wiw-l), Ai{t) = t"^ (40) 

to Eq. dnZI) with c^(e*^t) = \/ Ai{t) sin(6'/2) (since Pe (t) cx Ai{t) for a flxed background 
flow). As shown in Figure this corresponds to a constant shape, whose linear size grows 
hke with a 90° cusp at the rear stagnation point, where a branch point of ^^/w{w — 1) 
lies on the unit circle. For any flnite, Pe()f:), however, the cusp is smoothed, and the map 
remains univalent, although other singularities may form. Curiously, when mapped to the 
channel geometry with logz, the solution ()40|1 becomes a Saffman- Taylor flnger of width, 
A = 3/4. 
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3 Stochastic interfacial dynamics 



The continuous dynamics of conformal maps is a mature subject, but much attention is now 
focusing on analogous problems with discrete, stochastic dynamics. The essential change 
is in the kinematic condition: The expression for the interfacial velocity, e.g. Eq. is 
re-interpretted as the probability density (per unit arc length) for advancing the interface 
with a discrete "bump", e.g. to model a depositing particle. Continuous conformal-map 
dynamics is then replaced by rules for constru cting and composing the bum ps. This method 
of iterated conformal maps was introduced bv lHastings and Levitovl (|l998h in the context of 
Laplacian growth. 



Stoch astic Laplacian growth has been discussed since the early 1980s, but lHastings and Levitov 
( 19981 ) first showed how to implement it with conformal mapping. They proposed the fol- 
lowing family of bump functions, 



f\,e{w) 



fx {e-''w) , \w\>l 



w 



l-a 



'l + X){w + l] 
2w 



W+1 + \ w'^ + l-2w 



l + A 



(41) 
(42) 



as elementary univalent mappings of the exterior of the unit disk used to advance the interface 
(0 < a < 1). The function, fx^^w), places a bump of (approximate) area. A, on the unit 
circle, centered at angle, 6. Compared to analytic functions of the unit disk, the Hastings- 
Levitov function ()42|) generates a much more localized perturbation, foc used on the region 



betwe en two branch points, leaving the rest of the unit circle unaltered (jPavidovitch et al 
(1999)). For a = 1, the map produces a strike, which is a line segment of length 
emanating normal to the circle. For a = 1/2, the map is an invertible composition of simple 
linear, Mobius and Joukowski transformations, which inserts a semi-circular bump on the 
unit circle. As shown in Figure |H1 this yields a good descripti on of aggregating particle s, 
although other choices, like a = 2/3, have also been considered (jPavidovitch et al.l ()l999( )). 
Quantifying the effect of the bump shape remains a basic open question. 

Once the bump function is chosen, the conformal map, z = gn{w), from the exterior of the 
unit disk to the evolving domain with n bumps is constructed by iteration. 



9niw) = Qn-l if\„,e„iw)) 



(43) 



starting from the initial interface, given by go{w). All of the physics is contained in the 
sequence of bump parameters, {(A„,6'„)}, which can be generated in different ways (in the 
w plane) to model a variety of physical processes (in the z plane). As shown in Figure IH^b), 
the interface often develops a very complicated, fractal structure, which is given, quite re- 
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Figure 8: Simulation of the aggregation of (a) 4 and (b) 10,000 particles using the Hastings- 
Levitov algorithm (a = 1/2). Color contours show the quasi-steady concentration (or prob- 
ability) field for mobile particles arriving from infinity, and purple curves indicate lines of 
diffusive flux (or probability current). [Courtesy of Jaehyuk Choi and Benny Davidovitch.] 



markably, by an exact mapping of the unit circle. 



The great advantage of stochastic conformal mapping over atomistic or Monte Carlo simu- 
lation of interfacial growth lies in its mathematical insight. For example, given the sequence 
{(A„, 6n)} from a simulation of some physical growth process, the Laurent coefficients, Af:{n), 
of the conformal map , Qr,(w), as defined i n Eg. (|TH|l . can be calculated analytically. For the 
bump function (j42p . iDavidovitch et al.l (|l999[ ) provide a hierarchy of recursion relations, 
yielding formulae such as 



(44) 



m=l 



and explain how to interpret the Laurent coefficients. For example, Ai is the conformal radius 
of the cluster, a convenient measure of its linear extent. It is also the radius of a grounded 
disk with the same capacitance (with respect to infinity) as the cluster. The Koebe "1/4 
theorem" on univalent functions ( Duren ( IQSsi l) ensures that the cluster (image of the unit 
disk) is always contained in a disk of radius 4^41. The next Laurent coefficient, Aq, is the 
center of a uniformly charged disk, which would have the same asymptotic electric field as 
the cluster (if also charged). Similarly, higher Laurent coefficients encode higher multipole 
moments of the cluster. 
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Mapping the unit circle with a truncated Laurent expansion defines the web, which wraps 
around the growing tips and ex hibits a sort of turbulent dynamics, endlessl y forming and 
smoothing cusp- like protrusions (jHastingj (|l997l l: lHastings and LeyitoyI (|l998h ). The stochas- 
tic dynamics, however, does not suffer from finite-time singularities because the iterated map, 
by construction, remains univalent. In some sense, discreteness plays the role of surface ten- 
sion, as an another regularization of ill-posed continuum models like Laplacian growth. 

Diffusion-Limited Aggregation (PL A): T he stochastic analog of Laplacian growth is 
the DLA model of IWitten and Sander! ()l981l ). illustrated in Figure |H1 in which particles 
perform random walks one-by-one from inifinity until they stick irreversibly to a cluster, 
which grows from a seed at the origin. DLA and its variants (see below) provide simple 
models for many fractal patterns in nature, such as colloidal aggregates, dendritic electro- 
deposits, sn owfiakes, lightning s trikes , mineral deposits, and surface patterns in ion-beam 
microscopy (jBunde and Havlin In spite of decades of rese arch, however, DLA still 

presents theoretical mysteries, which are just beginning to unravel (jHalsevI (|200Cl )). 



The Hastings-Levitov algorithm for DLA prescribes the bump parameters, {(A„,^„)}, as 
follows. As in Laplacian growth, the harmonic function for the concentration (or probability 
density) of the next random walker approaching an n-particle cluster is simply. 



A Re log (7^^ ( 



(45) 



according to Eq. (jl)), since 4>{w) = ARe logw = log \w\ is the (trivial) solution to Laplace's 
equation in the mathematical w plane with = on the unit disk with a circularly symmetric 
fiux density. A, prescribed at infinity. Using the transformation rule, Eq. (j31|) . we then find 
that the evolving harmonic measure, Pn{z)\dz\, for the nth growth event corresponds to a 
uniform probability measure, Pn{9)d9, for angles. On, on the unit circle. 



w 



Pn{z)\dz\ = \V^(l)\\dz\ 



\9n-idw\ 



'Mdw] = — 
Zn 



Pn{0)d9, 



(46) 



where we set A = l/27r for normalization, which implicitly sets the time scale (see below). 
The conformal inv ariance of the harmonic mea sure is well known in mathematics, but the 
surprising result of iHastings and Levitovl (|l998( ) is that all the complexity of DLA is slaved 
to a sequence of independent, uniform random variables. 

Where the complexity resides is in the bump area, A„, which depends non-trivially on current 
cluster geometry and thus on the entire history of random angles, {6m\m, < n}. For DLA, 
the bump area in the mathematical w plane should be chosen such that it has a fixed value, 
Ao, in the physical z plane, equal to the aggregating particle area. As long as the new bump 
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is sufficiently small, it is natural to try to correct only for the Jacobian factor 



(47) 



m=l 



of the previous conformal map at the center of the new bump, 



An 



An 



■Jn-i(e 



(48) 



although it is not clear a priori that such a local approximation is valid. Note at least that 
(yf^ — >• oo, and thus A„ ^ 0, as the cluster grows, so this has a chance of working. 

Numerical simulations with the Hastings-Levitov algorithm do indeed produce nearly con- 
stant bump areas, as in Figure|Hl Nevertheless, much larger "particles", which fill deep fjords 
in the cluster, occasionally occur where the map varies too wildly, as shown in FigurelHfa). It 
is possible (but somewhat unsatisfying) to reject particles outside an "area acceptance win- 
dow" to produce rather realistic DLA clusters, as shown in Figure IHl^b). It seems that the 
rejected large b umps are so rare that they do not much influence statistical scaling properties 
of the clusters (|Stepanov and Levitovl (|200ll )). although this issue is by no means rigorously 
resolved. 

Fractal geometry: Fractal p atterns abound in nature, and DLA provides the most com- 
mon way to understand them iBunde and HavlinI ()l996[ ). The fractal scaling of DLA has 
been debated for decades, but conformal dynamics is shedding new light on the problem. 
Simulations show that the conformal radius (j44j) exhibits fractal scaling, Ai(n) oc n^^^f, 
where the fractal dimension, Dr = 1.71, agrees with the accepted value from Monte Carlo 



(random walk) simulations o: 



function (Davidovitch et al 



DLA, although the prefactor seems to depend on the bump 



igggr)). A perturbatlve renormalization-group analysis of the 



conformal dynamics by Hasting J ( 1997 ) gives a similar result, Df = 2 — 1/2-H/5 = 1 . 7. Th e 



multifractal spectrum of the harmonic measure has also been studied (|Jensen et al.l (|2002f l: 
Ball and Sor^Iil toO'Jh ). 



Perhaps the most basic question is whether DLA clusters are truly fractal - statistically self- 
similar and free of any length scale. This long-standing question requires accurate statistics 
and very large simulations, to erase the surprisingly long memory of the initial conditions. 
Conformal dynamics provides exact formulae for cluster moments, but simulations are limited 
to at most 10^ particles by poor O(n^) scaling, caused by the history-dependent Jacobian 
in Eq. (jlHl)- In contrast, efficient random- walk simulations can aggregate many millions of 
particles. 



Therefore. ISomfai et al.l (|l999| ) developed a hybrid method relying only upon the existence of 
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Figure 9: Simulations of fractal aggregates by IStepanov and Levitov ( 2001 ): (a) Super- 
imposed time series of the boundary, showing the aggregation of particles, represented by 
iterated conformal maps; (b) a larger simulation with a particle-area acceptance window; (c) 
the result of anisotropic growth probability with square symmetry; (d) square-anisotropic 
growth with noise control via flat particles; (e) triangular-anisotropic growth with noise 
control; (c) isotropic growth with noise control, which resembles radial viscous fingering. 
[Courtesy of Leonid Levitov.] 
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the conformal map, but not the Hastings-Levitov algorithm to construct it. Large clusters by 
Monte Carlo simulation, and approximate Laurent coefficients are computed, purely for their 
morphological information, as follows. For a given cluster of size N, M random walkers are 
launched from far away, and the positions, Zm, where they would ffist touch the cluster, are 
recorded. If the conformal map, z = (7„(e*^), were known, the points Zm would correspond 
to M angles 6^ on the unit circle. Since these must sample a uniform distribution, one 
assumes 9m = 21x171 /M for large M. From Eq. the Laurent coefficients are simply the 
Fourier coefficients of the discretely sampled function, Zm = J2 ^/cC*^™'^. Using this method, 
all Laurent coefficients appear to scale with the same fractal dimension. 



2/Df 



(49) 



although the ffist few coefficients crossover extremely slowly to the asymptotic scaling. 

Snowfiakes and viscous fingers: In conventional Monte Carlo simul ations, many vari- 
ants o f DLA have been proposed to model real patterns found in nature (jBunde and Havlin 
( 19961 )). For example, clusters closely ressembling snowfiakes can be grown by a combi- 
nation of noise control (requiring multiple hits before attachment) and anisotropy (on a 
lattice). Conformal dynamics offers the same flexibility, a s shown in Figure EI while a l- 
lowing anisotropy and noise to be controlled independently (|Stepanov and Levitovl (|2001l )). 
Anisotropy can be introduced in the growth probability with a weight factor, 1 + ccosma„, 
where a„ is the surface orientation angle in the physical plane given by Eq. (jl7|) . or by simply 
rejecting angles outside some tolerance from the desired symmetry directions. Noise can be 
controlled by flattening the aspect ratio of the bumps. Without anisotropy, this produces 
smooth fluid-like patterns (Figure IHl^f)), reminiscent of viscous fingers (Figure |21). 

The possible relation between DLA and viscous fingering is a tantalizing open question in 
pattern formation. Many authors have argued that the regularization of finite-time singular- 
ities in Laplacian growth by discreteness is somehow analogous to surface tension. Indeed, 
the average DLA cluster in a channel, grown by conformal m a pping , is similar (but not iden- 
tical) to a Saffman- Taylor finger of width 1/2 (jSomfai et al.l (j2003|)), and the instantaneous 
expected growth rate of a cluster can be related t o the Polubarinova-Galin (or "Shraiman- 
Bensimon" ) equation (jHastings and Levitovl (|l998h ) . Conformal dynamics with many bumps 
grown simultaneously suggests that Laplacian growth and DLA are in different universality 
classes, due to the b asic difference of layer-by-layer versus one-by-one growth, respectively 
( Barra et al.l (j2002a|)). Another multiple-bump algorithm with complete surface coverage, 
however, seems to yield the opposite conclusion (|Levermann and Procaccial (j200J)). 



Dielectric breakdown: In their original paper, iHastings and Levitovl (jl998|) allowed for 
the size of the bump in the physical plane to vary with an exponent, a, by replacing Jn-i with 



n-1 



.a/2 



in Eq. (gHD. In DLA 



a 



2), the bump size is roughly constant, but for < a < 2 
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Figure 10: Conformal-mapping simulations bv iHastingsl ( 2001 ) of the Dielectric Breakdown 
Model with (a) r/ = 2 and (b) r] = 3.5. [Courtesy of Matt Hastings.] 



the bump size grows with the local gradient of the Laplacian field. This is a simple model 
for dielectric breakdown, where the stochastic growth of an electric discharge penetrating a 
material is nonlinearly enhanced by the local electric field. One could use strikes (a = 0) 
rather than bumps (a = 1 /2) to better reproduce t he string- like branched patterns seen 
in laboratory experiments (|Bunde and HavlinI ()l996l )) and more familiar lightning strikes. 
The model displays a "stable-to-turbulent" phase transition: The relative surface roughness 
decreases with time for < a < 1 and grows for a > 1. 



The original Dielectric Breakdown Model (DBM) of iNiemever et al.l ()l984r ) has a more com- 
plicated conformal-dynamics representation. As usual, the growth is driven by the gradient 
of a harmonic function, 0, (the electrostatic potential) on an iso-potential surface (the dis- 
charge region). Unlike the a-model above, however, DBM growth events are assumed to 
have constant size, so the bump size in the mathematical plane is still chosen according to 
Eq. (|48|) . The difference lies in the growth measure, which does not obey Eq. (|46p . Instead, 
the generalized harmonic measure in the physical z plane is given by 



p{z) oc |V, 



(50) 



where t] is an exponent interpolating between the Eden model (?7 = 0), DLA (77 = 1), and 
nonlinear dielectric breakdown [rj > 1). For 77 7^ 1, the fortuitous cancellation in Eq. (jm]) 
does not occur. Instead, a similar calculation using Eq. yields a non-uniform probability 
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measure for the nth angle on the unit circle in the mathematical plane, 
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(51) 



which is complicated and depends on the entire history of the simulation. 



Nevertheless, conformal mapping can be applied fruitfully to DBM, because not solving 
Laplace's equation around the cluster outweighs the di fficultv of sampl ing the angle measure. 
Surmouting the latter with a Monte Carlo algorithm, iHastingsl (120011 ) has performed DBM 
simulations of 10"^ growth events, an order of magnitude beyond standard methods solving 
Laplace's equation on a lattice. The results, illustrated in Figure HUl support the theoretical 
conjecture that DBM clusters become one-dimensional, and thus non-fractal, for rj > 4. 

Using the conformal-mapping formalism, efforts are also underway to develop a unified scal- 
ing theory of the ?7-model for the growth probabi lity from DBM combined with the a-model 
above for the bump size (jBall and Somfail (|2002i l). 



Brittle fracture: Modeling the st ochastic dynamics of fract u re is a daunting problem, e s- 
pecially in heterogeneous materials (jBunde and HavlinI ()l996^ : iHermann and Roux jiggj). 
The basic equations and boundary conditions are still the subject of debate, and even the 
simplest models are difficult to solve. In two dimensions, stochastic conformal mapping 
provides an elegant, new alternative to discrete- lattice and finite-element models. 

In brittle fracture, the bulk material is assumed to obey Lame's equation of linear elasticity. 



(A + /i)V(V-u)+/iV2u 



(52) 



where u is the displacement field, p is the density, and p and A are Lame's constants. For 
conformal mapping, it is crucial to assume (i) two-dimensional symmetry of the fracture 
pattern and (ii) quasi-steady elasticity, which sets the left hand side to zero to obtain equa- 
tions of the type described above. For Mode III fracture, where a constant out-of-plane 
shear stress is applied at infinity, we have V ■ u = 0, so the steady Lame equation reduces 
to Laplace's equation for the out-of-plane displacement, V^Mx = 0, which allows the use of 
complex potentials. For Modes I and II, where a uniaxial, in-plane tensile stress is applied 
at infinity, the steady Lame equation must be solved. As discussed above, this is equivalent 
to the bi-harmonic equation for the Airy stress function, which allows the use of Goursat 
functions. 

For all three modes, the method of iterated conformal maps can be adapted to produce 



fracture patterns for a variety of physical assumptions about crack dynamics (jBarra et al 
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(2002b|)). For Modes I and II fracture, these models provide the first examples of stochastic 
bi-harmonic growth, which have interesting differences with stochastic Laplacian growth 
for Mode III fracture. The Hastings-Levitov formalism is used with constant-size bumps, 
as in DLA, to represent the fracture process zone, where elasticity does not apply. The 
growth measure a function of the excess tangential stress, beyond a critical yield stress, 
(Tc, characterizing the local strength of the material. Quenched disorder is easily included 
by making (Tc a random variable. In spite of its many assumptions, the method provides 
analytical insights, while obviating the need to solve Eq. ()52p during fracture dynamics, so 
it merits futher study. 

Advection-Diffusion-Limited Aggregation: Non-local fractal growth models typically 
involve a single bulk field driving the dynamics, such as the particle concentration in DLA, 
the electric field in DBM, or the strain field in brittle fracture, and as a result these models 
tend to yield statistically similar structures, apart from the effect of boundary conditions. 
Pattern formation in nature, however, is often fueled by multiple transport processes, such 
as diffusion, electromigration, and/or advection in a fluid flow. The effect of such dynam- 
ical competition on growth morphology is an open question, which would be difficult to 
address with lattice-based or flnite-element methods, since many large fractal clusters must 
be generated to fully explore the space and time dependence. 

Once again, conformal mapping provides a convenient means to formulate stochastic analogs 
of the non-Laplacian transport-limited growth models from §2.3 (in two dimensions). It is 
straightforward to adapt the Hastings-Levitov algorithm to construct stochastic dynamics 
driven by bulk flelds satsifying the conformally invariant system of equations ()35|) . A class 



of such models has recently been formulated bv lBazant et al.l (|2003^ . 



Perhaps the simplest case involving two transport processes, illustrated in Figure is 
Advection-Diffusion-Limited Aggregation (ADLA), or "DLA in a flow". Imagine a fluid car- 
rying a dilute concentration of sticky particles flowing past a sticky object, which begins to 
collect a fractal aggregate. As the cluster grows, it causes the fluid to flow around it and 
changes the concentration fleld, which in turn alters the growth probability measure. Assum- 
ing a quasi-steady potential flow with a uniform speed far from the cluster, the dimensionless 
transport problem is 

Pe oV0 ■ Vc = V^c, = 0, z G n,{t) (53) 

c = 0, n-V0 = O, a = n-Vc, z e dn,{t) (54) 
c — s> 1, V0 X, \z\ oo (55) 

where Peo is the initial Peclet number and a is the diffusive flux to the surface, which 
drives the growth. The transport problem is solved in the mathematical w plane, where it 
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Figure 11: A simulation of Advection-Diffusion-Limited Aggregation from iBazant et al 
( 2003[ ) ■ In each row, the growth probabihties in the physical z-plane (on the right) are 
obtained by solving advection-diffusion in a potential flow past an absorbing cylinder in the 
mathematical w-plane (on the left), with the same time-dependent Pelcet number. 
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corresponds to a uniform potential flow of concentrated fluid past an absorbing circular cylin- 
der. The normal diffusive flux on the cylinder, a{6,Pe), ca n be obtained from a tabulated 



numerical solution or an accurate analytical approximation (|Choi et al.l (j2004b|)). 



Because the boundary condition on at infinity is not conformally invariant, the flow in 
the w plane has a time-dependent Peclet number, Pe (t) = Ai{t)PeQ, which grows with the 
conformal radius of the cluster. As a result, the probability of the nth growth event is given 
by a time-dependent, non-uniform measure for the angle on the unit circle, 

Pn(^n) = I- r„ a(e*^", (56) 

where /3 is a constant setting the mean growth rate. The waiting time between growth events 
is an exponential random variable with mean, r^, given by the current integrated flux to the 
object, 

^= r a{e\A,{t^_,))de. (57) 

P^n Jo 

Unlike DLA, the aggregation speeds up as the cluster grows, due to a larger cross section to 
catch new particles in the flow. 

As shown in Figure ^2 f he model displays a universal dynamical crossover from DLA (the 
unstable fixed point) to an advection-dominated stable fixed point, since Pe (t) oo. Re- 
markably, the fractal dimension remains constant during the transition, equal to the value 
for DLA, in spite of dramatic changes in the growth rate and morphology (as indicated by 
higher Laurent coefficients). Moreover, the shape of the "average" ADLA cluster in the high- 
Pe regime of Figure [TT] is quite similar (but not identical) to the exact solution, Eq. (j^ . 
for the analogous continuous problem in Figure [7| Much remains to be done to understand 
these kinds of models and apply them to materials problems. 



4 Curved Surfaces 



Entov and Etingoa (|44| ) considered the generalized problem of Hele-Shaw flows in a non- 
planar cell having non-zero curvature. In such problems, the velocity of the viscous flow 
is still the (surface) gradient of a potential, (p, but this function is now a solution of the 
so-called Laplace- Beltrami equation on the curved surface. The Riemann mapping theorem 
extends to curved surfaces and says that any simply-connected smooth surface is conformally 
equivalent to the unit disk, the complex plane, or the Riemann sphere. A common example 
is the well-known stereographic projection of the surface of a sphere to the (compactified) 
complex plane. Under a conformal mapping, solutions of the Laplace-Beltrami equation map 
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Figure 12: Conformal- mapping simulation of DLA on a sphere. Particles diffuse one by one 
from the North Pole and aggegrate on a seed at a South Pole. [Courtesy of Jaehyuk Choi, 
Martin Bazant, and Darren Crowdy.] 



to solutions to Laplace's equation and this combination of facts led lEntov and Etingoa (j4 



to identify classes of explicit solutions to t he continuous Hele- Shaw problem in a variety of 



non-planar cells. With very similar intent. iParisio et al.l (|2001[ ) have recently considered the 



evolution of Saffman- Taylor fingers on the surface of a sphere. 

By now, the reader may realize that most of the methods already considered in this article 
are, in principle, amenable to generalization to curved surfaces, which can be reached by 
conformal mapping of the plane. For example. Figure ^1 shows a simulation of a DLA 
cluster growing on the surface of a sphere, using a generalized Hastings-Levitov algorithm, 
which takes surface curvature into account. The key modification is to multiply the Jacobian 
in Eq. (P7|) by the Jacobian of the stereographic projection, 1 + where R is the radius 

of the sphere. 

It should also be clear that any continuous or discrete growth model driven by a conformally- 
invariant bulk field, such as ADLA, can be simulated on general curved surfaces by means 
of appropriate conformal projection to a complex plane. The reason is that the system 
of equations ()35|) is invariant under any conformal mapping, to a flat or curved surface, 
because each term transforms like the Laplacian, J'V'^cp, where J is the Jacobian. 
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The purpose of studying these models is not only to understand growth on a particular ideal 
shape, such as a sphere, but more generally to explore the effect of local surface curvature on 
pattern formation. For example, this could help interpret mineral deposit patterns in rough 
geological fracture surfaces, which form by the diffusion and advection of oxygen in slwoly 
flowing water. 



5 Outlook 

Although conformal mapping has been with us for centuries, new developments with ap- 
plications continue to the present day. This appears to be the first pedagogical review of 
stochastic conformal-mapping methods for interfacial dynamics, which also covers the latest 
progress in continuum methods. Hopefully, this will encourage the futher exchange of ideas 
(and people) between the two fields. Our focus has also been on materials problems, which 
provide many opportunites to apply and extend conformal mapping. 

Building on specific open questions scattered throughout the text, we close with a general 
outlook on directions for future research. A basic question for both stochastic and continuum 
methods is the effect of geometrical constraints, such as walls or curved surfaces, on interfacial 
dynamics. Most work to date has been for either radial or channel geometries, but it would 
be interesting to describe finite viscous fingers or DLA clusters growing near walls of various 
shapes, as is often the case in materials applications. 

The extension of conformal-map dynamics to multiply connected domains is another math- 
ematically challenging area, which has received some attention recently but seems ripe for 
further development. Understanding the exact solution structure of Laplacian-growth prob- 
lems using the mathematical abstraction of quadrature domain theory holds great potential, 
especially given that mathematicians have already begun to e xplore the exten t to which the 



various mathematical concepts extend to higher- dimensions (jShapird ()l992f l). Describing 
multiply connected domains could pave the way for new mathematical theories of evolving 
material microstructures. Topology is the main difference between an isolated bubble and 
a dense sintering compact. Microstructural evolution in elastic solids may be an even more 
interesting, and challenging, direction for conformal-mapping methods. 

From a mathematical point of view, much remains to be done to place stochastic conformal- 
mapping methods for interfacial dynamics on more rigorous ground. This has recently been 
achieved in the simpler case of Stochastic Loe wner Evolution (SLE) , which has a similar his- 



tory to the interfacial problems discussed here (|Kager and Nienhuid (j200J)). Oded Schramm 



introduced SLE in 2000 as a stochastic version of the continuous Loewner evolution from 
univalent function theory, which grows a one- dimensional random filament from a disk or 
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half plane. This important development in pure mathematics came a few years after the pio- 
neering DLA papers of Hastings and Levitov in physics. A notable difference is that SLE has 
a rigorous mathematical theory based on stochastic calculus, which has enabled new proofs 
on the properties of percolation clusters and self-avoiding random walks (in two dimensions, 
of course). One hopes that someday DLA, DBM, ADLA, and other fractal-growth models 
will also be placed on such a rigorous footing. 

Returning to materials applications, it seems there are many new problems to be considered 
using conformal mapping. Relatively little work has been done so far on void electromigra- 
tion, viscous sintering, solid pore evolution, brittle fracture, electrodeposition, and solifica- 
tion in fluid flows. The reader is encouraged to explore these and other problems using a 
powerful mathematical tool, which deserves more attention in materials science. 
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